home *** CD-ROM | disk | FTP | other *** search
/ APDL Eductation Resources / APDL Eductation Resources.iso / programs / astronomy / skyview / !SkyView / c / Riset < prev    next >
Encoding:
Text File  |  1993-08-26  |  4.7 KB  |  139 lines
  1. /*------------------------------------------------------*/
  2. /*                         Riset                        */
  3. /*    Set of functions for use in the calculation of    */
  4. /*   rising, setting and culmination times of objects   */
  5. /*          with time-dependent RA and Dec.             */
  6. /*------------------------------------------------------*/
  7.  
  8. #include <math.h>
  9.  
  10. #include "os.h"
  11. #include "menu.h"
  12. #include "coords.h"
  13. #include "sv_header.h"
  14. #include "riset.h"
  15. #include "datime.h"
  16.  
  17. /*------------------------------------------------------*/
  18. /*   Set maximum declination, beyond which object is    */
  19. /*   considered so close to the celestial pole that it  */
  20. /*                   does not move.                     */
  21. /*------------------------------------------------------*/
  22. #define RISET_NOMOVE (REAL)1.57
  23.  
  24.  
  25. /*------------------------------------------------------*/
  26. /* Function which attempts to find time of culmination, */
  27. /*               given an initial guess.                */
  28. /*------------------------------------------------------*/
  29.   BOOL riset_cul(observerstr *ob_ptr, int maxit, int converge,
  30.                  radecfntype *radecfn, int id, int *hour2ptr, int *min2ptr)
  31. {
  32.   int hour1, min1;  /* Time1.                           */
  33.   int diff;         /* Difference between time1 & time2.*/
  34.   REAL ra, dec;     /* RA & Dec of object.              */
  35.   int iteration=0;  /* Iteration No.                    */
  36.  
  37. /* Iterate to find time of culmination.                 */
  38.   do
  39.   {
  40.     iteration +=1;
  41.     hour1 = *hour2ptr;
  42.     min1  = *min2ptr;
  43.   /* Calculate RA and Dec of object at time1.           */
  44.     radecfn(id, ob_ptr, hour1, min1, &ra, &dec);
  45.   /* Get time2, when this bit of sky culminates.        */
  46.     datime_time_today(ra, ob_ptr, hour2ptr, min2ptr);
  47.     diff = 60 * (*hour2ptr - hour1) + *min2ptr - min1;
  48.     if (diff < 0) diff = -diff;
  49.   }
  50.   while (iteration < maxit  &&  diff > converge);
  51.  
  52. /* Inform caller whether or not the time converged.     */
  53.   return (diff <= converge);
  54. }
  55.  
  56.  
  57. /*------------------------------------------------------*/
  58. /*   Function which attempts to find time of rising or  */
  59. /*          setting, given an initial guess.            */
  60. /*------------------------------------------------------*/
  61.   BOOL riset_riset(observerstr *ob_ptr, int maxit, int converge,
  62.                    radecfntype *radecfn, int id, BOOL rise, REAL horalt,
  63.                    int *hour2ptr, int *min2ptr)
  64. {
  65.   int hour1, min1;  /* Time1.                           */
  66.   int diff;         /* Difference between time1 & time2.*/
  67.   REAL ra, dec;     /* RA & Dec of object.              */
  68.   int iteration=0;  /* Iteration No.                    */
  69.   double d_dec;     /* double version of declination.   */
  70.   double coshrang;  /* Cos of hour angle at 'horizon'.  */
  71.   double hourang;   /* Hour angle at 'horizon'.         */
  72.   BOOL cos_valid;   /* Is coshrang a valid cosine?      */
  73.   REAL lst;         /* Local siderial time @ rise/set.  */
  74.   /* double version of observer's latitude:             */
  75.   double d_latit = (double)ob_ptr->latit;
  76.  
  77.  
  78. /* Iterate to find time of rising / setting.            */
  79.  
  80.   do
  81.   {
  82.     iteration +=1;
  83.     hour1 = *hour2ptr;
  84.     min1  = *min2ptr;
  85.  
  86.   /* Calculate RA and Dec of object at time1.           */
  87.     radecfn(id, ob_ptr, hour1, min1, &ra, &dec);
  88.  
  89.   /* Give up if declination is too close to 90 degrees. */
  90.     if (dec >= RISET_NOMOVE || dec <=-RISET_NOMOVE)  return FALSE;
  91.  
  92.   /* Get cosine of hour angle when patch of sky         */
  93.   /* specified by ra and dec rises.                     */
  94.     d_dec = (double)dec;
  95.     coshrang = (sin(horalt) - sin(d_latit)*sin(d_dec)) /
  96.                (cos(d_latit)*cos(d_dec));
  97.  
  98.   /* Convert cosine to angle.  Note that for rising,    */
  99.   /* true hour angle is actually 2PI - hourang.         */
  100.     if      (coshrang >  1.0)
  101.     {
  102.       if (iteration == maxit) return FALSE;
  103.       hourang = 0.0;
  104.       cos_valid = FALSE;
  105.     }
  106.     else if (coshrang < -1.0)
  107.     {
  108.       if (iteration == maxit) return FALSE;
  109.       hourang = DblPI;
  110.       cos_valid = FALSE;
  111.     }
  112.     else
  113.     {
  114.       hourang = acos(coshrang);
  115.       cos_valid = TRUE;
  116.     }
  117.  
  118.   /* Get local siderial time when patch of sky has this */
  119.   /* hour angle.                                        */
  120.     if (rise)
  121.       lst = ra - (REAL)hourang;
  122.     else
  123.       lst = ra + (REAL)hourang;
  124.  
  125.   /* Get time2, corresponding to this lst.              */
  126.     datime_time_today(lst, ob_ptr, hour2ptr, min2ptr);
  127.  
  128.   /* Calculate difference between time2 & time1.        */
  129.     diff = 60 * (*hour2ptr - hour1) + *min2ptr - min1;
  130.     if (diff < 0) diff = -diff;
  131.  
  132.   }
  133.   while ( (diff > converge || !cos_valid) && iteration < maxit );
  134.  
  135. /* Inform caller if time did not converge or cosine of  */
  136. /* hour angle was invalid.                              */
  137.   return (diff <= converge && cos_valid);
  138. }
  139.